
# Replication file for Appendix G.1.: Figure 12, Tables 9 and 10


#install.packages("Hmisc")
require(Hmisc)
#install.packages("stargazer")
require(stargazer)

load("Brexit_Rep.RData")

leavers.sub <- brexit.sub2[which(brexit.sub2$leaver == 1),]
remainers.sub <- brexit.sub2[which(brexit.sub2$leaver == 0),]

dvs <- Cs(migs.take.jobs.6,migs.more.terror.6,closed.immigrants.6,hurt.standing.refugee.6,
          threaten.culture.refugee.6,overwhelm.welfare.refugee.6)

# ATEs for Leavers (create empty vector to populate)
ate.l <- vector(length(dvs), mode = "list")
names(ate.l) <- dvs

# ATEs for Remainers
ate.r <- vector(length(dvs), mode = "list")
names(ate.r) <- dvs

for (i in 1:(length(dvs))){
  
  modelformula <- paste(dvs[i],"~ post")
  print(modelformula)
  
  ate.l[[dvs[i]]] <- eval(substitute(lm(.modelformula, data = leavers.sub, 
                                                             weights = w8w6), list(.modelformula = modelformula)))
  
}


for (i in 1:(length(dvs))){
  
  modelformula <- paste(dvs[i],"~ post")
  print(modelformula)
  
  ate.r[[dvs[i]]] <- eval(substitute(lm(.modelformula, data = remainers.sub, 
                                             weights = w8w6), list(.modelformula = modelformula)))
  
}



cis.l <- matrix(NA,length(dvs),2)
b.l <- matrix(NA,length(dvs),1)
cis.r <- matrix(NA,length(dvs),2)
b.r <- matrix(NA,length(dvs),1)


for (i in 1:length(dvs))
{
  b.l[i] <- ate.l[[i]]$coefficients[2]
  ci <- confint(ate.l[[i]])
  lb <- ci[2,1]
  ub <- ci[2,2]
  cis.l[i,] <- cbind(lb,ub)
}

for (i in 1:length(dvs))
{
  b.r[i] <- ate.r[[i]]$coefficients[2]
  ci <- confint(ate.r[[i]])
  print(b.r[i])
  lb <- ci[2,1]
  ub <- ci[2,2]
  cis.r[i,] <- cbind(lb,ub)
}


y.lab <- seq(from = 1, to=length(dvs),by=1)
vnames <-c(migs.take.jobs.6 = "Migs. Take\nJobs",migs.more.terror.6 = "Migs. Bring\nTerror", closed.immigrants.6 = "Not\nOpen to\nMigs.",
           hurt.standing.refugee.6 = "Refs. Don't\nImprove\nUK Image", threaten.culture.refugee.6 = "Refs.\nThreaten\nCulture", overwhelm.welfare.refugee.6 ="Refs.\nOverwhelm\nServices")


jitter <-  0.1

#pdf("Fig1_NoControls.pdf", height = 5, width=5)
par(mar =c(3,3.1, 4.1, 0.01), oma=c(1,5,.2,.2), cex.axis=.8, cex.main=.9)
plot(b.l, y.lab, type="n", ylab ="", xlab = "", yaxt="n", xlim = c(-0.5,0))
axis(2, at=y.lab, labels=vnames, las = 2, cex.axis=0.8)
mtext("Average Treatment Effects",side=1,line=2,outer=F, cex = 0.8)
points(b.l, y.lab, pch= 19, cex=.8, col = "black")
segments(cis.l[,1], y.lab, cis.l[,2], y.lab, lty=1, col = "black")
segments(cis.r[,1], y.lab+jitter, cis.r[,2], y.lab+jitter, lty=1, col = "black")
points(b.r, y.lab+jitter, pch= 19, cex=.8, col = "white")
points(b.r, y.lab+jitter, pch= 1, cex=.8)
abline(v=0)
legend(-0.5, 6, pch=c(19,1), legend=c("Leave", "Remain"), cex =.8, border = "black")
dev.off()

### Tables

stargazer(ate.l[[1]], ate.r[[1]], ate.l[[2]], ate.r[[2]], ate.l[[3]], ate.r[[3]])
stargazer(ate.l[[5]], ate.r[[5]], ate.l[[6]], ate.r[[6]], ate.l[[4]], ate.r[[4]])








